{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Introduction to BIDMat"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "BIDMat is a multi-platform matrix library similar to Matlab, Julia or Numpy/Scipy. Its intended primarily for machine learning, but is has a broad set of operations and datatypes and should be suitable for many other applications. \n",
    "\n",
    "BIDMat is probably unique in its integration of CPU and GPU data types. Other features include:\n",
    "* Interactivity. Thanks to the Scala language, BIDMat is interactive and scriptable. \n",
    "* Performance, thanks to CPU and GPU native code, and to Scala's speed on memory-bound operations. \n",
    "* Parallelism, thanks to Scala's actor framework and parallel collection classes. \n",
    "* Rich, open syntax of math operators, +,-,*,/,⊗,∙,∘\n",
    "* Runs on JVM, extremely portable, leverages Hadoop, Yarn, Spark etc. \n",
    "\n",
    "BIDMat is a library that is loaded by a startup script, and a set of imports that include the default classes and functions. We include them explicitly in this notebook. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 CUDA device found, CUDA version 7.0\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(0.84340703,10187849728,12079398912)"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import BIDMat.{CMat,CSMat,DMat,Dict,IDict,FMat,GMat,GIMat,GSMat,GSDMat,HMat,IMat,Image,LMat,Mat,SMat,SBMat,SDMat}\n",
    "import BIDMat.MatFunctions._\n",
    "import BIDMat.SciFunctions._\n",
    "import BIDMat.Solvers._\n",
    "import BIDMat.JPlotting._\n",
    "\n",
    "Mat.checkMKL\n",
    "Mat.checkCUDA\n",
    "Mat.setInline\n",
    "if (Mat.hasCUDA > 0) GPUmem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These calls check that CPU and GPU native libs loaded correctly, and what GPUs are accessible.\n",
    "If you have a GPU and CUDA installed, GPUmem will printout the fraction of free memory, the absolute free memory and the total memory for the default GPU. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Basic Matrix Algebra"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From this cell onward, we are in the BIDMat environment. Let define some matrices and basic algebra on them. BIDMat has Float, Double, Int and Complex matrix types. We'll start with integer matrices. To create an array of ones, do"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1   1   1   1\n",
       "   1   1   1   1\n",
       "   1   1   1   1\n",
       "   1   1   1   1\n"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val a = ones(4,4)   // ones creates a 4x4 Float matrix (FMat). FMat is the default.\n",
    "                    // you use prefixes to get other types. e.g. iones gives an integer matrix. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BIDMat.FMat\n"
     ]
    }
   ],
   "source": [
    "%type a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can create row and column matrices (FMat) by listing their elements:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "2,2,4,4,5,5"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "row(2,2,4,4,5,5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   5\n",
       "   4\n",
       "   3\n",
       "   2\n"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "col(5,4,3,2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also create integer row or column matrices (IMat) with irow and icol. These functions accept lists of values, *or* ranges:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "0,1,2,3,4,5,6,7,8,9"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val r = irow(0 until 10)   // until gives a range excluding the last element"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BIDMat.IMat\n"
     ]
    }
   ],
   "source": [
    "%type r"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0\n",
       "   1\n",
       "   2\n"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "icol(0 to 2)               // to gives a range including the last element"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can create a matrix of sequential integers like this"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4   8  12\n",
       "   1   5   9  13\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val b = izeros(4,4)             // An integer matrix this time, filled with zeros\n",
    "b(?) = icol(0 until 16)         // Now treat as a 1d array, fill with a range\n",
    "b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The questionmark ? is BIDMat's wildcard character. Even though b is two-dimensional, b(?) linearizes its contents into a 16-element column and puts the RHS into it. \n",
    "\n",
    "The RHS should be another 16x1 integer matrix (IMat), but when supplied with a range (0 until 16), BIDMat automatically casts the range to an IMat. This is called an implicit conversion in Scala.\n",
    "\n",
    "From the order of elments in the array after the assignment, you can see that BIDMat uses *Column-major-order*. This is similar to Matlab, Fortran and Julia, but different from C and Python which are row-major. \n",
    "\n",
    "Transpose is implemented with a \"t\" operator:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   1   2   3\n",
       "   4   5   6   7\n",
       "   8   9  10  11\n",
       "  12  13  14  15\n"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val bt = b.t"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Basic Math Operators"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Math operators have their expected results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1   5   9  13\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n",
       "   4   8  12  16\n"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val c = a + b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BIDMat.FMat\n"
     ]
    }
   ],
   "source": [
    "%type c             // We added an integer matrix (IMat) and a float matrix (FMat), so what type is the result?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "BIDMat implicitly casts IMats to FMats to perform algebraic operations. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  -1   3   7  11\n",
       "   0   4   8  12\n",
       "   1   5   9  13\n",
       "   2   6  10  14\n"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b - a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   6  22  38  54\n",
       "   6  22  38  54\n",
       "   6  22  38  54\n",
       "   6  22  38  54\n"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a * b                    // Matrix multiply (not element-wise multiply)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4   8  12\n",
       "   1   5   9  13\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b / a                    // This is element-wise division, some toolkits instead multiply by the inverse of a. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Advanced Math Operators"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As well as these operators, BIDMat includes several other important operators with their standard unicode representation. They have an ASCII alias in case unicode input is difficult. Here they are:\n",
    "\n",
    "<pre>\n",
    "Unicode operator    ASCII alias    Operation\n",
    "================    ===========    =========\n",
    "       ∘                *@         Element-wise (Hadamard) product\n",
    "       ∙                dot        Column-wise dot product\n",
    "       ∙→              dotr        Row-wise dot product\n",
    "       ⊗               kron        Kronecker (Cartesian) product\n",
    "</pre>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4   8  12\n",
       "   1   5   9  13\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b ∘ a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "6,22,38,54"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b ∙ a "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  24\n",
       "  28\n",
       "  32\n",
       "  36\n"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b ∙→ a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   0   0   0   4   4   4   4   8   8   8   8  12  12  12  12\n",
       "   0   0   0   0   4   4   4   4   8   8   8   8  12  12  12  12\n",
       "   0   0   0   0   4   4   4   4   8   8   8   8  12  12  12  12\n",
       "   0   0   0   0   4   4   4   4   8   8   8   8  12  12  12  12\n",
       "   1   1   1   1   5   5   5   5   9   9   9   9  13  13  13  13\n",
       "   1   1   1   1   5   5   5   5   9   9   9   9  13  13  13  13\n",
       "   1   1   1   1   5   5   5   5   9   9   9   9  13  13  13  13\n",
       "   1   1   1   1   5   5   5   5   9   9   9   9  13  13  13  13\n",
       "  ..  ..  ..  ..  ..  ..  ..  ..  ..  ..  ..  ..  ..  ..  ..  ..\n"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b ⊗ a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> TODO: using the operators above, construct a 5x5 matrix such that every element is one greater than the element to the left, and the element above. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Transposed Multiplies"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Matrix multiply is the most expensive step in many calculations, and often involves transposed matrices. To speed up those calcualtions, we expose two operators that combine the transpose and multiply operations:\n",
    "\n",
    "<pre>\n",
    "^&ast;  - transpose the first argument, so a ^&ast; b is equivalent to a.t &ast; b\n",
    "&ast;^  - transpose the second argument, so a &ast;^ b is equivalent to a &ast; b.t\n",
    "</pre>\n",
    "these operators are implemented natively, i.e. they do not actually perform transposes, but implement the effective calculation. This is particulary important for sparse matrices since transpose would involve an index sort. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   6  22  38  54\n",
       "   6  22  38  54\n",
       "   6  22  38  54\n",
       "   6  22  38  54\n"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a ^* b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   6  22  38  54\n",
       "   6  22  38  54\n",
       "   6  22  38  54\n",
       "   6  22  38  54\n"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a.t * b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  24  24  24  24\n",
       "  28  28  28  28\n",
       "  32  32  32  32\n",
       "  36  36  36  36\n"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b *^ a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  24  24  24  24\n",
       "  28  28  28  28\n",
       "  32  32  32  32\n",
       "  36  36  36  36\n"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b * a.t"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Edge Operators"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Most operators support scalar arguments as in the last section. There are also many situations where its helpful to apply an operation with an \"edge\" argument, that is a vector whose long dimension matches the matrix. This is similar to the \"broadcast\" feature in Numpy. For example, we can define a vector v as:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1\n",
       "   2\n",
       "   3\n",
       "   4\n"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val v = col(1,2,3,4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1   1   1   1\n",
       "   2   2   2   2\n",
       "   3   3   3   3\n",
       "   4   4   4   4\n"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v ∘ a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The elementwise multiply by v is applied to every column of a. We could also apply v.t to the rows of a:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1   2   3   4\n",
       "   1   2   3   4\n",
       "   1   2   3   4\n",
       "   1   2   3   4\n"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v.t ∘ a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Edge operators arise in several context, e.g. normalizing a matrix along rows or columns. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The sum() function computes the sums of columns of b and returns a row vector:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "10,26,42,58"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sum(c)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can use this to normalize b along its columns:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  0.10000  0.19231  0.21429  0.22414\n",
       "  0.20000  0.23077  0.23810  0.24138\n",
       "  0.30000  0.26923  0.26190  0.25862\n",
       "  0.40000  0.30769  0.28571  0.27586\n"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val d = c / sum(c)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "1,1,1,1"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sum(d)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "or rows"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  0.035714   0.17857   0.32143   0.46429\n",
       "  0.062500   0.18750   0.31250   0.43750\n",
       "  0.083333   0.19444   0.30556   0.41667\n",
       "   0.10000   0.20000   0.30000   0.40000\n"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val e = c / sum(c,2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1\n",
       "   1\n",
       "   1\n",
       "   1\n"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sum(e,2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    ">TODO: Using the sqrt() function on matrices, normalize the columns of the matrix c so that their L2 norm (or equivalently the dot product of the column with itself) is 1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "ename": "",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "incomplete"
     ]
    }
   ],
   "source": [
    "//val cn =  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "ename": "",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "incomplete"
     ]
    }
   ],
   "source": [
    "//cn ∙ cn "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Scalars"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Scalars (primitive numerical values) can be used in most expressions:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   2   2   2   2\n",
       "   2   2   2   2\n",
       "   2   2   2   2\n",
       "   2   2   2   2\n"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a + 1                      // Add an integer, result is still an FMat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  3.5000  3.5000  3.5000  3.5000\n",
       "  3.5000  3.5000  3.5000  3.5000\n",
       "  3.5000  3.5000  3.5000  3.5000\n",
       "  3.5000  3.5000  3.5000  3.5000\n"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val aa = 3.5 * a             // Multiply by a double, result is still an FMat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BIDMat.FMat\n"
     ]
    }
   ],
   "source": [
    "%type aa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Note** In Scala, floating point numbers like 3.5 have double precision by default. Many languages would cast the last result to a double matrix since double is the smallest container for both Float and Double data. We argue that when someone writes 3.5 x a, they mean to scale the matrix a by that factor and preserve its type, not to cause a type conversion. Single-precision constants in Scala need an \"f\" suffix, i.e."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  3.7000  3.7000  3.7000  3.7000\n",
       "  3.7000  3.7000  3.7000  3.7000\n",
       "  3.7000  3.7000  3.7000  3.7000\n",
       "  3.7000  3.7000  3.7000  3.7000\n"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a + 2.7f"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In case you encounter double matrices in a calculation without meaning to, it may because of operations with double-precision constants. Use the floating point notation above for scalars to minimize the chances of unintentional up-conversion to double matrices.\n",
    "\n",
    "For the next step, we will need a floating point version of the matrix b, which we can construct like this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4   8  12\n",
       "   1   5   9  13\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val c = FMat(b)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    ">TODO: Create a float vector of values from -10 to 10 spaced by 0.1. Then apply the logistic function 1/(1+exp(-c)) to and call the plot() function on the results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "ename": "",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "<console>:2: error: illegal start of simple expression\r",
      "       val y=\r",
      "       ^"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    }
   ],
   "source": [
    "val x=\n",
    "val y=\n",
    "plot(x,y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Indexing and Slicing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You access invidual array elements using parenthesis (unlike [] in Python)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1   1   1   1\n",
       "   1   2   1   1\n",
       "   1   1   1   1\n",
       "   1   1   1   1\n"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a(1,1) = 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1   1   1   1\n",
       "   1   2   1   1\n",
       "   1   1   1   1\n",
       "   1   1   1   1\n"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "1.0"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a(0,3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can use the wildcard ? to access rows and columns:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4   8  12\n",
       "   1   5   9  13\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   4\n",
       "   5\n",
       "   6\n",
       "   7\n"
      ]
     },
     "execution_count": 45,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(?,1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "2,6,10,14"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(2,?)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ranges work as expected:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   5   9\n",
       "   6  10\n"
      ]
     },
     "execution_count": 47,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(1->3, 1->3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And you can use arbitrary integer vectors to access submatrices:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4  12\n",
       "   1   5  13\n"
      ]
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(icol(0,1), icol(0,1,3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Another shorthand constructor for integer matrices is the backslash operator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "0,1,3"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val ii = 0\\1\\3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BIDMat.IMat\n"
     ]
    }
   ],
   "source": [
    "%type ii"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And this syntax is handy for indexing expressions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4  12\n",
       "   1   5  13\n"
      ]
     },
     "execution_count": 51,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(0\\1, 0\\1\\3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Slices can be used for assignment:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   4   8  12\n",
       "   0   0   9   0\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n"
      ]
     },
     "execution_count": 52,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(1, 0\\1\\3) = 0\\0\\0\n",
    "b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and you can use scalars on the RHS to simplify bulk assignments"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  -1  -1   8  -1\n",
       "  -1  -1   9  -1\n",
       "   2   6  10  14\n",
       "   3   7  11  15\n"
      ]
     },
     "execution_count": 53,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(0\\1, 0\\1\\3) = -1\n",
    "b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Matrices also accept single indices that reference elements in column-major order:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "7"
      ]
     },
     "execution_count": 54,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(7)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also use vectors of indices or ranges to assign arbitrary elements of a matrix, or all of them:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "16,15,14,13,12,11,10,9,8,7,6,5,4,3,2,1"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b(0->16) = (16 to 1 by -1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> TODO: Define a set of indices ii such that for any 4x4 matrix m, m(ii) = m.t. \n",
    "\n",
    "> **HINT:** you already computed it!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0.62465  0.050579   0.49750   0.52764\n",
       "   0.63891   0.73763   0.59721  0.059120\n",
       "   0.84862   0.39953   0.80420   0.95167\n",
       "   0.26990   0.91040   0.55897   0.51923\n"
      ]
     },
     "execution_count": 56,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val m = rand(4,4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "ename": "",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "incomplete"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    }
   ],
   "source": [
    "// val ii = "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "ename": "",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "incomplete"
     ]
    }
   ],
   "source": [
    "// m(ii)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Reducers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "BIDMach has several \"reducers\" that aggregate along rows or columns. We already saw one of these, which was the <code>sum()</code> function. Two other important ones are maxi() and mini(). These both compute the max or min respectively, along columns (default) or rows. e.g."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0.086997    0.10430    0.13608    0.39235   0.032846\n",
       "    0.58714    0.20012    0.93311    0.94488    0.59928\n",
       "    0.21865    0.23329    0.30459  0.0090942    0.60729\n",
       "    0.41737    0.67991    0.71646    0.11169   0.067964\n"
      ]
     },
     "execution_count": 59,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val x = rand(4,5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "0.58714,0.67991,0.93311,0.94488,0.60729"
      ]
     },
     "execution_count": 60,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val xm = maxi(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "  0.39235\n",
       "  0.94488\n",
       "  0.60729\n",
       "  0.71646\n"
      ]
     },
     "execution_count": 61,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val xm2 = maxi(x,2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Its often very useful to know not only what element was the max or min, but also its *index*. The functions <code>maxi2</code> and <code>mini2</code> do this. They have a \"2\" suffix to indicate that they return 2 values. The first is the max or min *value*, the second is the max or min *index*:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "0.58714,0.67991,0.93311,0.94488,0.60729"
      ]
     },
     "execution_count": 62,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val (vmax, imax) = maxi2(x)\n",
    "vmax"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "1,3,1,1,2"
      ]
     },
     "execution_count": 63,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "imax"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The first 3 means that the max element of the first column was in row number 3, etc. We can similarly compute the min along rows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0.032846\n",
       "    0.20012\n",
       "  0.0090942\n",
       "   0.067964\n"
      ]
     },
     "execution_count": 64,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val (vmin, imin) = mini2(x,2)\n",
    "vmin"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   4\n",
       "   1\n",
       "   3\n",
       "   4\n"
      ]
     },
     "execution_count": 65,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "imin"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A last important reducer is <code>accum</code> which is similar to Matlab's <code>accumarray</code>, or numpy's <code>accum</code>. It allows you to tally some values into specific positions in an output array. The format is:\n",
    "\n",
    "<pre>\n",
    "accum(inds, vals, nrows, ncols)\n",
    "</pre>\n",
    "\n",
    "where inds is an nx2 matrix (IMat) of row,column indices, vals are the values to sum there, and nrows and ncols are the matrix dimensions. Its easiest to see with an example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0   0\n",
       "   0   1\n",
       "   2   3\n",
       "   0   1\n",
       "   3   3\n"
      ]
     },
     "execution_count": 66,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val inds = 0\\0 on 0\\1 on 2\\3 on 0\\1 on 3\\3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 67,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1\n",
       "   2\n",
       "   3\n",
       "   4\n",
       "   5\n"
      ]
     },
     "execution_count": 67,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val vals = col(1f, 2f, 3f, 4f, 5f)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   1   6   0   0   0\n",
       "   0   0   0   0   0\n",
       "   0   0   0   3   0\n",
       "   0   0   0   5   0\n"
      ]
     },
     "execution_count": 68,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "accum(inds, vals, 4, 5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can see that each of the <code>vals</code> was saved in the position specified by the corresponding <code>inds</code>. Most of the locations to save were distinct, except for the second and fourth rows, which specified the same location. Those two values (2f and 4f) were summed in that location. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Find Functions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The <code>find</code> function is similar to Matlab's find, and Numpy's <code>nonzero</code> function. It comes in several flavors, depending on how many values are returned:\n",
    "\n",
    "<pre>\n",
    "val ii = find(m)            // find and return the single indices of non-zero elements of m\n",
    "val (ii, jj) = find2(m)     // find and return the (row, column) indices of non-zeros in m\n",
    "val (ii, jj, vv) = find3(m) // find an return the (row, column) indices and values of non-zeros of m\n",
    "</pre>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    ">TODO: For the matrix below, use find3 to deconstruct it into row, column and value matrices. Then use accum to build it up again. You can use the \\ operator to horizontally concatenate two matrices."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "   0.70932   0.68385   0.88133   0.48188   0.85245\n",
       "   0.53094   0.33278   0.24578   0.66883   0.50532\n",
       "   0.68560   0.14570  0.077930   0.10479   0.74823\n",
       "   0.78596   0.82603   0.67407   0.36298   0.86334\n"
      ]
     },
     "execution_count": 69,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val rr = rand(4,5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "incomplete"
     ]
    }
   ],
   "source": [
    "// deconstruct rr\n",
    "// val rr2 =                              // Rebuild it"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## For and While Loops"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We used ranges before. There are two flavors, closed or \"to\" ranges, and open or \"until\" ranges."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "Range(0, 1, 2, 3, 4, 5)"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "0 to 5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "Range(0, 1, 2, 3, 4)"
      ]
     },
     "execution_count": 72,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "0 until 5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For loops use ranges in a natural way"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "run number 0\n",
      "run number 1\n",
      "run number 2\n",
      "run number 3\n",
      "run number 4\n"
     ]
    }
   ],
   "source": [
    "for (i <- 0 until 5) {\n",
    "    println(\"run number %d\" format i)\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "While loops provide no special loop variable management, simply a test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "counting down 6\n",
      "counting down 5\n",
      "counting down 4\n",
      "counting down 3\n",
      "counting down 2\n",
      "counting down 1\n"
     ]
    }
   ],
   "source": [
    "var i = 6\n",
    "while (i > 0) {\n",
    "    println(\"counting down %d\" format i)\n",
    "    i -= 1\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Warning: Performance Sinkhole!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For loops are much more complex that while loops. They apparently create a local evaluation context for each iteration, and the overhead is several times higher than for while loops. Lets measure this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 75,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(-39.68925,2.729)"
      ]
     },
     "execution_count": 75,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import java.util.Random\n",
    "val randgen = new Random\n",
    "tic\n",
    "var sum = 0f\n",
    "for (i <- 0 until 1000000) {\n",
    "    sum += randgen.nextFloat - 0.5f\n",
    "}\n",
    "val t = toc\n",
    "(sum, t)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 76,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(-361.9767,1.644)"
      ]
     },
     "execution_count": 76,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "tic\n",
    "var sum = 0f\n",
    "var i = 0\n",
    "while (i < 1000000*100) {\n",
    "    sum += randgen.nextFloat - 0.5f\n",
    "    i += 1\n",
    "}\n",
    "val t = toc\n",
    "(sum, t)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For loops (and Scala's other functional tools like sequence classes) are extremely powerful for e.g. multi-threading and multi-GPU computing. But they're not suitable for lightweight iteration over elements. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Functional Programming in BIDMat"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Functional programming in BIDMat (or Numpy or Matlab) avoids explicit iteration over the elements of matrices and concentrates instead on whole-array operations and (if irregular access is needed) on manipulation of index matrices. Its not unlike the use of global operations on DataFrames. \n",
    "\n",
    "This approach allows highly-parallelized code to be used to implement these routines. It often makes for more succinct code and (with some practice), greater readability. We'll concentrate on applying those ideas in the next part of the Lab."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## I Feel the Need for Speed !!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Yes, its not only fun to have a fast toolkit, but it really matters for performance. Not just runtime, but most algorithms can trade off time for precision by simply training more thoroughly, or training a richer model. \n",
    "\n",
    "BIDMat/BIDMach is currently the only system which fully integrates GPU computing with CPU computing. Its only the only system to have fully rooflined sparse matrix primitives. This is very important, since these are the bottleneck for machine learning on the most common types of data (text, web, clickthrough etc). Let's measure exactly how much difference this makes.\n",
    "\n",
    "First, we'll define a few matrices, both on the CPU and on a GPU (assuming you have one)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(     0,     0)     2\n",
       "(    40,     0)     1\n",
       "(    69,     0)     1\n",
       "(   177,     0)     1\n",
       "(   679,     0)     1\n",
       "(   715,     0)     1\n",
       "(  1802,     0)     1\n",
       "(  3368,     0)     1\n",
       "    ...    ...    ...\n"
      ]
     },
     "execution_count": 77,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val n = 8192\n",
    "val a = rand(n,n)                     // a random dense matrix (CPU)\n",
    "val b = powrand(n,n,100)              // a random power-law sparse matrix with 100 elements per column (CPU)\n",
    "val ga = GMat(a)                      // a GPU version of a\n",
    "val gb = GSMat(b)                     // a GPU version of b\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we could just go ahead and do our calculations on a,b,ga,gb directly. This is a common scenario. But we would also like to illustrate BIDMat/BIDMach's support for generics. So instead will create variables of type \"Mat\" to hold those variables, and perform arithmetic on those instead."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "null"
      ]
     },
     "execution_count": 78,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var ma:Mat = a              // create a generic Mat variable and bind it to a\n",
    "var mb:Mat = b              // create a generic Mat variable and bind it to b\n",
    "var mc:Mat = null           // we'll use this to hold results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we'll benchmark both dense and sparse matrix multiply. Dense CPU multiply first."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(97.716995,11.252)"
      ]
     },
     "execution_count": 79,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "flip\n",
    "mc = ma * ma\n",
    "gflop"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The \"flip\" function starts a timer and reset the flop count. \"gflop\" returns two values: the gigaflop count, and the time since the last \"flip\". "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(3.637876,2.348)"
      ]
     },
     "execution_count": 80,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "flip\n",
    "mc = ma * mb\n",
    "gflop"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now lets bind those variables to GPU matrices instead:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(     0,     0)     2\n",
       "(    40,     0)     1\n",
       "(    69,     0)     1\n",
       "(   177,     0)     1\n",
       "(   679,     0)     1\n",
       "(   715,     0)     1\n",
       "(  1802,     0)     1\n",
       "(  3368,     0)     1\n",
       "    ...    ...    ...\n"
      ]
     },
     "execution_count": 81,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ma = ga\n",
    "mb = gb"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and run exactly the same code:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 82,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(3088.5159,0.356)"
      ]
     },
     "execution_count": 82,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "flip\n",
    "mc = ma * ma\n",
    "gflop"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You'll probably see a good order-of-magnitude speedup over the CPU calculation. This shouldnt be surprising. GPUs have a well-earned reputation for dense-matrix performance. \n",
    "\n",
    "What's less well-known is their sparse matrix performance, which yields roughly order-of-magnitude gains as well."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 83,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(40.67492,0.21)"
      ]
     },
     "execution_count": 83,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "flip\n",
    "mc = ma * mb\n",
    "gflop"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should see performance in the 20-50 gflops range, which is near the roofline for sparse operations on the current generation of Nvidia GPUs. \n",
    "\n",
    "This is very important, because that operation (and two other variants that have similar performance) is the dominant step in most common machine learning algorithms. With careful design of the entire learning pipeline, you can translate that advantage into a end-to-end speedup by the same factor. \n",
    "\n",
    "Furthermore, by writing generic code (using the Mat class as above) you can hide the details of implementation (CPU vs GPU) and run your algorithm in either environment. You can also support either sparse or dense matrices, and many of BIDMach's learning algorithms will work with either."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> TODO: try the transpose operator &ast;^ in the cells above. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Matrix Caching"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One of the challenges of working with GPUs is the current lack of memory management (i.e. a garbage collector). The very high streaming speed of GPU memory makes it very difficult to do memory management without significant slowdowns. BIDMach includes a matrix caching scheme which allows re-use of matrix storage. It works particularly well in BIDMach's minibatch algorithms, which process same-sized blocks of data many times. \n",
    "\n",
    "To understand caching, lets first notice that every matrix has a unique long id or guid:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 84,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(40341906236551069,4295235833250888165)"
      ]
     },
     "execution_count": 84,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val a = rand(4,4)\n",
    "val b = rand(4,4)\n",
    "(a.GUID, b.GUID)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Normally, when you do calculations with matrices, new containers are created to hold the results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 85,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(-5272250351793888220,-1073144993423657928)"
      ]
     },
     "execution_count": 85,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val c = a * b\n",
    "val d = a + b\n",
    "(c.GUID, d.GUID)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "But with caching enabled, the same expression will yield the same container (different expressions with the same arguments will still yield different containers):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 86,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "(-3486012619384765122,-3486012619384765122,4427540394353907121)"
      ]
     },
     "execution_count": 86,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Mat.useCache = true\n",
    "val c = a + b\n",
    "val d = a + b\n",
    "val e = a * b\n",
    "(c.GUID, d.GUID, e.GUID)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Although this approach causes aliasing (c and d now point to the same container), in a functional programming language the same expressions should always hold the same value. Arrays are mutable objects, so are not guaranteed to hold the same value. Nevertheless if you program in functional style, the same expressions should hold the same value and caching is a safe operation. \n",
    "\n",
    "With that caveat, caching is a very helpful performance optimization. The learners in BIDMach automatically turn caching on and off when you run a learning algorithm and in this way are able to eliminate matrix allocation after the first iteration. Its a necessary feature to be able to use a GPU on large datasets, and it often accelerates calculations on the CPU by removing memory allocation and garbage collection overhead. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 87,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": []
    },
    {
     "data": {
      "text/plain": [
       "false"
      ]
     },
     "execution_count": 87,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Mat.useCache = false"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "name": "scala",
   "version": "2.11.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
